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Abstract 



Motivated by recent Monte-Carlo simulations of Hoglund and Sandvik (arXiv:0808.0408), we 
study edge response in square lattice quantum antiferromagnets. We use the 0(3) non-linear 
cr-model to compute the decay asymptotics of the staggered magnetization, energy density and 
local magnetic susceptibility away from the edge. We find that the total edge susceptibility is 
negative and diverges logarithmically as the temperature T — > 0. We confirm the predictions of 
the continuum theory by performing a 1/5 expansion of the microscopic Heisenberg model with 
the edge. We propose a qualitative explanation of the edge dimerization seen in Monte-Carlo 
simulations by a theory of valence-bond-solid correlations in the Neel state. We also discuss the 
extension of the latter theory to the response of a single non-magnetic impurity, and its connection 
to the theory of the deconfined critical point. 
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I. INTRODUCTION 



The Heisenberg antiferromagnet on a square lattice is one of the best known model 
magnetic systems. It has been studied extensively both numerically by quantum Monte- 
Carlo and analytically by 1/S expansion and field-theoretic methods. It is known to have 
an ordered ground state at zero temperature with the staggered magnetization reduced by 
quantum fluctuations to N b = (N) = 0.307 for the spin 5* = 1/2. 1 

Despite many years of study, the simple Heisenberg model does not cease to surprise us. 
Recent Monte-Carlo simulations 2 on the S = 1/2 model have shown that the edge response 
in this system is very peculiar. In particular, a negative edge susceptibility is observed at 
low temperatures. This result is in contrast with an intuitive picture of a "dangling" edge 
spin leading to an enhancement in the susceptibility. The simulation of local susceptibility 
near the edge shows that the negative sign does not come from the edge spins per se, whose 
susceptibility is, indeed, enhanced, but rather from a tail in the response decaying away from 
the edge. Another curious effect observed in Ref. [2] is the dimerization of bond response near 
the edge, leading to the appearance of a comb-like structure, as in Fig. [T] The tendency to 
dimerize into singlets near the edge was argued in Ref. [2] to be the source of negative edge 
susceptibility. 



FIG. 1: A schematic picture of the comb structure in bond strengths observed in Monte-Carlo 
simulations^, with a free edge on the left side. 



In the present paper, we study large-distance asymptotics of the edge response of a 
square lattice quantum antiferromagnet by means of an effective 0(3) a- model description. 
This field-theoretic method is an expansion in powers of energy and momentum, with the 
microscopic physics entering at each order through a finite number of parameters, such as 
the spin- wave speed c, the spin stiffness p s and the value of the staggered moment N b . 1 
The 0(3) a-model has proved powerful for studying finite temperature/size effects, which 
typically lead to a crossover into an 0(3) model of lower dimensionP It turns out to be also 



We will use the subscript b from here on to denote bulk properties. 
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useful for studying the edge behaviour, particularly as no new parameters beyond the bulk 
ones are needed to describe the leading low temperature, large distance asymptotics in the 
edge response. We concentrate our attention on the staggered moment (N(x)), the local 
energy density (e(x)) and the local magnetic susceptibility X-L^)- We show that at zero 
temperature these quantities approach their bulk values away from the edge with simple 
power law forms, 



(N(x)) - N b 



Nb 8irp s x 



X±(x)~X±,b = —z (1.3) 

Sttxc 



Integrating eq. (1.3), we conclude that the total edge susceptibility per unit edge length is 



negative and diverges logarithmically with the system size 

1 



X±,ed ge = -g-^ log(L/°) C 1 - 4 ) 



We show that at finite temperature the 1/x power law in the susceptibility (1.3) is cut- 
off for distances larger than the thermal wave-length, x > c/T, leading to the total edge 
susceptibility, 

X±,ed gc = q — log(c/Ta) (1.5) 
one 

Such a log divergent susceptibility is indeed seen in the Monte Carlo simulations^. For the 
co-efficient of the logarithm in Xedge = (2/3)x_i_,edge, with c = 1.69 J, we find —0.0157/ J, 
while the Monte Carlo has a best fit value of —0.0182/ J (see Fig. [2|. This is in reasonable 
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FIG. 2: Edge susceptibility: Comparison of the Monte Carlo data of Ref. |2] (dots) with the best fit 
line Jxedgc = -0.0182 log(0. 219 J/T) to the low T data. 
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agreement, with the difference probably attributable to difficulties in numerically reaching 
the asymptotic low T limit. 

As for the edge comb structure seen in Ref. |2| this is a short distance phenomenon, 
which cannot be studied within our continuum 0(3) a-model. In fact, the standard, "per- 
turbative" treatment of the 0(3) model describes only the low-energy excitations which live 
near the wave-vector (tt, tt) and cannot provide any information about valence-bond-solid 
correlations, which live near (tt, 0) and (0,7r). Because these correlations are gapped in the 
antiferromagnet, they must decay exponentially away from the edge, as seen in Monte-Carlo. 
To capture the short-distance physics, we have performed a 1/5* expansion of the Heisenberg 
model on the lattice with an edge. We find the large-distance asymptotics in agreement with 
the predictions of our continuum theory. However, we don't reproduce the multiple short- 
distance oscillations of bond energies away from the edge seen by Monte-Carlo. Instead, we 
find that the bonds touching the edge are stronger than the bulk, while all the subsequent 
bonds are weaker. We conclude that the edge dimerization is, likely, a non-perturbative 
effect in 1/S, which is invisible in the spin- wave expansion. It is remarkable that such non- 
perturbative effects are present in the simple 5 = 1/2 Heisenberg model, where the 1/S* 
expansion yields quantitatively accurate results for many quantities. 

In principle, one may be able to explicitly incorporate the non-perturbative physics in 
the form of hedgehogs into the semi-classical, large 5* treatment of the Heisenberg model. 
The hedgehog configurations are relevant for the dimerization physics, as they carry Berry 
phases,^ which endow them with non-trivial quantum numbers under the lattice symmetryP 
However, studying the hedgehog contribution to the edge physics is technically intractable. 

Instead, we pursue a more phenomenological approach, in which we assume that the 
system possesses a dynamical valence-bond-solid order parameter with a large correlation 
length. This assumption is justified close to a phase transition into a valence-bond-solid 
phase, which can be tuned by adding additional frustrating interactions to the Heisenberg 
modelP^ Moreover, even for the pure, nearest neighbour Heisenberg model with S = 1/2, it 
has been argued long agd^ that the quantum fluctuations are strong enough that the system 
is "proximate" to a phase transition at which the magnetic order is lost. This proximity 
is manifested by the existence of an intermediate temperature window, dominated by the 
quantum critical point (the low temperature physics is dominated by the antiferromagnet, 
while the high temperature physics is dominated by the non-universal lattice effects). The 
observation of edge dimerization over more than 5 lattice spacings in the latest Monte Carlo 
simulations implies that the correlation length of the valence-bond-solid order parameter in 
the S = 1/2 Heisenberg model is rather large, further supporting the proximity to a phase 
transition. 

We show that the comb structure of the bond order seen in Monte-Carlo simulations 
can be qualitatively understood in the quantum critical language. The particular details 
of the critical theory are not very important for this purpose - the physics can be read 
off straight-forwardly from the transformation properties of observables under the lattice 
symmetry. In particular, we demonstrate that close to the critical point the oscillations of 
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bonds perpendicular to the edge and lines parallel to the edge in the comb can be related 
to each other. 

In another recent paper with Kaul and Melko^, we have discussed the response of the 
valence bond solid order to a single non-magnetic impurity (such as a Zn site replacing a Cu 
site). We used there a phenomenological theory similar in spirit to that used here for the edge 
response. We will review that theory here and also discuss its connection to the impurity 
response in the deconfmed theory of the Neel to valence bond solid transition discussed in 
Ref. HOl For this single-site impurity case, we are able to infer the non-perturbative role of 
hedgehogs and Berry phases in somewhat greater detail. 

This paper is organized as follows. Section II A is devoted to the description of the edge 
in the framework of the 0(3) model at zero temperature. In section II B we discuss the 
crossover of edge susceptibility to finite temperature. In section III we perform the large S 
expansion of the Heisenberg model with an edge. In section IV we discuss edge dimerization 
in a quantum antiferromagnet in the proximity to a phase transition into a valence-bond- 
solid. Finally, in section [V] we discuss the related theory of the response near a non-magnetic 



impurity. Some concluding remarks are presented in section VI 



II. EDGE RESPONSE IN THE 0(3) cr-MODEL 
A. Zero Temperature 

In this section we discuss the large distance asymptotic behaviour away from the edge of 
the staggered moment, local uniform susceptibility and the bond energies using the contin- 
uum 0(3) cr-model. The advantage of this approach is that the results obtained are exact, 
depending only on a few phenomenological parameters, such as spin-wave velocity c, spin- 
stiffness p s and bulk staggered moment Nj,. These parameters are known from l/S'-expansion 
and Monte-Carlo simulations. 

The a-model action for the local order parameter ft, satisfying n 2 = 1, is 

S=j J d 2 xdr(d^n) 2 (2.1) 

Here, we've set c = 1, we will restore c at the end of the computations. Since we are studying 
the problem with an edge, we also have to consider boundary perturbations. The simplest 
terms allowed by symmetries are, 

Abound = X) W dydr (d^n) 2 (2.2) 

This term is irrelevant by power counting (the coupling has scaling dimension —1), and can 
be ignored for the leading asymptotic behaviour calculations performed below. Note that 
the "lower dimension" surface term nd x n vanishes identically due to the constraint n 2 = 1. 
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The absence of a boundary term, implies that n obeys free boundary conditions, 



d x n = (2.3) 



as can be seen by varying the action (2.1) with respect to n, integrating by parts and 
requiring that the surface term be zero. 

To set up perturbation theory, we write n = (n, y/l — tt 2 ) and expand the action in jf, 
obtaining, 

S = & J d 2 xdr {{d^f + ^l^^d^) 2 ) (2.4) 

The second term in brackets above can be expanded as a power series in n - yielding terms 
with couplings of scaling dimension —1 and lower. These terms again will not influence the 
leading asymptotic behaviour of observables discussed below. 

We are, thus, left with the free theory for the Goldstone fields n, supplemented by the 
free boundary condition d x n = 0. The propagator with these boundary conditions is, 

(-K a (x,r)A^,r')) = S 4 I 2 T I T r2 e^-T') e ^iy-y) cos{kxX ) cos (^) 

p u s J In In 7r u z + k l x + 

5 ab 

= —(D(x - x',y - y',r - t>) + D(x + x',y - y',r - t')) (2.5) 

Ps 

where D(x) is the standard 3d massless propagator, 

D(x) = (2.6) 
4:7r\x\ 

Now, we can calculate the observables. Let's start with the staggered moment (N). 
The microscopic N(x) is related to the 0(3) field n(x) via a multiplicative renormalization, 

— # 

N{x) = Nj,Z]yn(x) where Nf, is the exact value of the bulk staggered magnetization and Zjy 
is a formal power series in p" 1 , adjusted order by order to give (iV 3 ) = N b in the bulk. 
Hence, the staggered moment, to leading order is, 

(n 3 (x)} = (1 - f ) = 1 - 1(£>(0) + D(2x, 0, 0)) = 1 - 1(£>(0) + ±-) (2.7) 

Thus, as lim^oo Z^{m?{x)) = 1, and to leading order p° = p s , 

1 . 1 f d 3 k 1 



ZK = 1 + - m = 1 + -J—- (2.8) 
which is the familiar expression known from calculations with no boundary. So, 



"^l 1 -^,! (2 - 9) 
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FIG. 3: Depletion of the staggered moment, —((N(x)) — iVj,) near the edge. The dotted line is 
the calculation in the 1/S expansion. The solid line is the 0(3) a-model result for asymptotic 
behaviour, with phenomenological parameters p s , c, N b matched to 1/5" expansion. 



where we've reinserted the spin-wave velocity c. The result (2.9) is asymptotically exact 



and shows suppression of the Neel moment near the edge. We can check the result (2.9) 



against the large distance asymptotics of the 1/S* expansion performed in section III The 
parameters p s , c and Nb are known in 1/S expansion to be at leading order, 



Ps = JS 2 



2V2JSa, N b = S 



(2.10) 



where a is the lattice spacing. Substituting these parameters into (2.9) and comparing to 



our numeric integration results from 1/S" expansion on the lattice with an edge, we find very 
reasonable asymptotic agreement (see Fig. [3]). 

Next we consider the uniform transverse susceptibility x±- Recall, the uniform magnetic 



field enters (2.1 ) as, 



S 



H 



d 2 xdr ((d T n a - ie abc H b n c f + (c^n) 2 ) 



(2.11) 



The corresponding response function is, 

X a \x,x>) = 5H f { ^ b{xf) = P° s (S ab - (n a n\x)))5\x - x') - {p G s fe acd e^{n c d T n d {x)n^{x')) 

(2.12) 
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Specializing to the transverse susceptibility, a, b = 1, 2 and expanding in tF, 

X ab (x,x') « ^(^-(^(^(xO))^^-^)^-^) 

- (p°)V c e 6d ((<9 r 7r c (x)<9V(z')) + ((<9 r 7T c (x)(7T d 7f<9 r 7f - V^Tr^x')) + (x <-> x', c <-> d))) 

(2.13) 

Now, we are actually interested in local response to a static, uniform external field, 

xf (x) = lim ! dx'dy'dT'x ah {x,x')e- iqy{v ~ v,) (2.14) 

q y ^0 J 

Note that for a finite system size/temperature relevant for Monte-Carlo simulations, at zero 
external field, there is no distinction between parallel and transverse susceptibility, and we 
expect, 

X(x) = \x±{x) (2.15) 
Since we are working with the static susceptibility, the contribution of the terms in the last 



two lines of (2.13) is zero, and 



X f(x) = p° s (6 ab - (n a (x)n\x'))) = P y\l - \(D(0) + L>(2x,0,0))) (2.16) 

Ps 



We know that in the bulk, x±,b — ^ m x^ooX±i x ) = Ps by Lorentz invariance. The bare 
spin-stiffness p° s = p s Z p where Z p is a formal power series in l/p s . Thus, 



1 . , 1 f d 3 k 



Z„ = l + -D(0) = l + - I — (2.17) 



and we recognize the standard renormalization factor for p s . Note that the equality of the 
first non-trivial terms in and Z p is an accident, which occurs in the 0(3) model (for 
0(N) the coefficients are generally different). Thus, 

Xx{x) = ^-^- (2.18) 

where we've reinserted c. Note that the deviation of X_i_(x) from its bulk value is negative, 
in agreement with Sandvik's simulations. Moreover, the long distance contribution to the 
total edge susceptibility (per edge length) is given by, 

f°° 1 
X±,ed ge = / dx(x±(x) - x±,b) —\og(L x /a) (2.19) 

J0 87TC 

At zero temperature, the log divergence of the long-distance tail will always overpower any 
short-distance contribution (which can be positive as suggested by the 1/5 calculation in 



section III), leading to a negative total edge susceptibility, as seen by Sandvik. At a finite 
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temperature T (and in the infinite volume limit) the log L x divergence will be cut-off at the 
"thermal length," cT -1 , leading to 

Xx,*~-^log(^) (2.20) 

This result will be confirmed by an explicit calculation in the next section. 

Finally, we come to the behaviour of the bond energies. We observe that the sum of 
bonds energies along the x and y directions is just the local energy density 

e(x) ~ -^(SiS i+& + SiS i+ y) (2.21) 
For the free field theory describing our Goldstones, in Minkowski space, 

<*) = | ((^) 2 + (W 2 ) (2-22) 
Continuing this to Euclidean space, 

<*) = f ("(^) 2 + (^) 2 ) (2-23) 

Now, 

-^-(d tl 7f(x)d u 7f(x)) = Hm q lu ( D ( x ~x',y- y', T-T')+D(x + x',y- y', r-r')) (2.24) 

The first term on the righthandside is independent of the distance from the edge and, 
therefore, we drop it. Noting, 

^W = -J^(^-3^) (2.25) 



the second term in (2.24) yields 



((<9 r vff (x)> = -d 2 D(2x,0,0) = - * .„ (2.26) 



2 u-r , v-// ~r v 47r ( 2x )3 

|((^vf) 2 (x)) = +^(2^,0,0) = ^-y 3 (2.27) 

J <(M a (^)> = -d 2 y D(2x,0,0) = (2.28) 

Collecting terms we obtain 



(<*)) = J7T— 3 (2-29) 
Note that energy density is enhanced near the edge, corresponding to a decrease of bond 



strengths, —(SiSj). We can again compare the asymptotically exact expression (2.29) to 
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FIG. 4: Asymptotic increase of local bond energy near the edge. The dotted line is the calculation 
in the 1/5" expansion. The solid line is the 0(3) cr-model result for asymptotic behaviour, with 
phenomenological parameters p s , c, iV& matched to 1/5 expansion. 



the results of the 1/S expansion in section III by using the parameters (2.10). We see from 
Fig. [4] that the agreement is rather good. 



B. Edge susceptibility at finite temperature 

To compute the uniform susceptibility at finite temperature T < /) s , we follow the usual 
strategy of dividing the field n(x, r) into zero frequency piece, n(x) and finite frequency 
modes n a (x, r), 

(2.30) 



n a (x , r) = y/1 - n a n a n a (x) + n a (x, r)e a a {x) 

where a = 1, 2 and e a (x) and n(x) form an orthonormal basis. The strategy is to first inte- 
grate over the "fast" modes ir a to obtain an effective action for the slow n field. Expanding 
the action in powers of 7r to leading order, 

S « p l J d 2 xdr (d^ a ) 2 + ^ J d 2 xdr {{d in a f{\ - n 2 ) + d^d^n? + 2die%e^ a d^ p ) 

(2.31) 

In setting up the perturbation theory in 7r the first term above is treated as the free piece, 
while the coupling of tt to the slow fields in the second term is treated as a perturbation. 
Thus, in a theory with the edge at finite temperature, the bare propagator for the ir field 
still satisfies free boundary conditions, 



(7r Q (f,r)7T /3 (f / ,r / )) = —8 a pD n {x,x') 

Ps 



(2.32) 
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where, 

D n (x, x') = D(x -x',y- y', r - r) + D(x +x',y- y', r - r') (2.33) 



with 

72 



(2-34) 



Now, expanding the susceptibility (2.12), 



X a \x) = p° s (5 ab - (n a n b (x))) - (p° s ) 2 e acd e b ^ / dV^^eftf K^x)^^')) 

(2.35) 



At leading order, we may factorize the correlator of slow e and fast tc fields in (2.35). 
Moreover, since at finite temperature rotational invariance is restored, 

Sjab Sjab 

^n b {x)) = °—(n 2 {x)) = °— (2.36) 

Hence, the local susceptibility becomes, 

X ab (x) = P° s 5 ab - (p° s ) 2 e acd e be f J dV(e^(x)e^(f' ^d^x^d^x')) (2.37) 

We see that the susceptibility involves a convolution of correlators of slow and fast fields. 
Evaluating the correlation function of the fast fields explicitly, 

x ab (x) = 2 -p^-e acd e^{8 ai 8p S -8 aS 8 Pl ) J d 3 x\e c a e$(x)e e 7 e f 5 (x'))(d T D n (x, x')f 

(2.38) 

We note, 

J dr'(d T D n (x,x')) 2 = ^ ^^£> n (f,f ,u n ) 2 
= — ^2 u n(^(^ ~ % \ u n ) 2 + 2D(x - x, uo n )D(x - Rx, u n ) + D(x - Raf , oo n ) 2 ) 

(2.39) 

where R denotes reflection across the edge at x — 0. In the absence of an edge, we can 



drop the last two terms in (2.39). Then we note that the correlation function of it's decays 



exponentially for large distances, hence only \x — x'\ < T 1 contribute to the integral in 



(2.37). The slow degrees of freedom n(x) and e a (x) fluctuate only on much larger distances 
(in fact T _1 serves as an effective short-distance cut-off for the slow degrees of freedom), 
hence we can to leading order set x = x' in the correlation function of the e's. This leads to 
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a considerable simplification as, 



e a a e b a = 5 ab - n a n° 



and, 



and 



X a \x) 



5 n />5 



aSVp-y 



)(e c a e d (x)ete f s (x)) = \{Sr°Sfi - 5^5 de ) 



5 ab ( P° 



d 6 x\d T D n (x,x')) 2 



(2.40) 



(2.41) 



(2.42) 



Now let's introduce the edge back. We wish to compute the deviation of local susceptibility 



from its bulk value. The major difference from the situation in the bulk is that eq. (2.39) 



no longer depends just on the difference x — x'. For xT < 1, the integral over x! in (2.38) is 
saturated with x'T < 1 and hence, we can effectively set x — x' = 0, y = y' in the correlation 



function of the e's and recover the simple form (2.42). However, for xT 3> 1, the part of the 



integral in (2.38) that represents x( x )~Xb is no longer saturated at x' ~ x. Hence, one really 
has to compute the correlation function of the slow degrees of freedom. For T _1 
we expect this to modify x( x ) ~ Xb (which, as we shall see, is exponentially suppressed as 
e -4irTx^ logarithmic corrections. On the other hand, for x > £, we expect additional 
exponential suppression coming from the slow degrees of freedom. As we shall see, the total 



edge susceptibility is saturated by xT < 1 and, hence, can be computed directly from (2.42). 



Keeping the above remarks in mind, we obtain from (2.39) and (2.42) 



X{ x ) 



P° s - 2^ £ u£ / dx' / dy'{D{ 

P u, n ^0 J ~°° J -°° 



(2.43) 



The first term under the integral in (2.43) is the familiar temperature dependent correction 



to bulk susceptibility, while the second term represents the edge contribution. Performing 
the integral over x', 



4 1 



X (x) = Xb(T) - -- 



UJ 2 e 2ik x x 



3 ^ (2^)2(^+^2)2 



(2.44) 



where, 



Xb(T) 



d 2 k uo 2 n 



2 P. 

3 c 2 2np s 



(2.45) 



Now, we can compute the asymptotics of (2.44). For xT/c <C 1, we can replace the sum 
over ui n by an integral, 



X{x) -> Xb{T) 



4 f d 3 k 



„ 2 2 & 



3 J (2vr) 3 (k 2 + uj 2 f 



Xb{T) - 



i r d 2 k i 



j 2ft A>;"r ~3^j 



3 J (27T) 2 k 



Xb(T) - 



Ylnxc 
(2.46) 
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which agrees with our earlier T = result (2.18) upon the usual replacement (2.15). In the 



opposite limit xT/c ^> 1, the sum in (2.44) is going to be dominated by the smallest thermal 

mass, ui n= i, and, 

, , 2T 
3^ 



Xb 



xT\ 2 



A-kTx/c 



(2.47) 



As noted earlier, this result will be modified by logarithmic corrections for x <C £ and 



additional exponential suppression for x £. It is also now clear from (2.47) that the total 



edge susceptibility is saturated by xT < 1, so that the corrections mentioned above can be 



ignored for its computation, and we can use eq. (2.44), which obeys the scaling form, 



XO) ~Xb = T fx( Tx ) 



(2.48) 



Thus, 



X 



dx(x{x) - Xb) 



duf(u) 



(2.49) 



Ta 



where a is a short distance cut-off. We observe that the singular behaviour of Xedge for T — > 
can be extracted from the short distance asymptotic of x( x ) (2.46). Noting, f x (u) 1 
for u 



Yl-KU 



0. 



Xcdg 



du 



12tt J Ta u 

as predicted from T = behaviour in the previous section. 



— logf-) 
12ttc & \TaJ 



(2.50) 



III. LARGE S EXPANSION OF THE HEISENBERG MODEL WITH AN EDGE 

In this section we perform the large S expansion of the Heisenberg model on a square 
lattice with an edge. We start with the usual nearest neighbour Hamiltonian, 

H = jJ2SiSj (3.1) 

(ij) 

and use the Holstein-Primakoff representation of spin operators, which at leading order in 
1/5 reads, 



3f = S- b\bi, St = V2Sb t , Sr = V2Sb\, % e A (3.2) 
SI = -S + c\a, St = V2Sc\, s; = V2Sc h % e B (3.3) 

where A and B are the two sublattices. We place the edge at i x = 0. Utilizing the transla- 
tional invariance along the y direction, 

I ^ 
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where —tt/2 < k y < n/2 and N y is the number of sites in the y direction, we obtain the 
Hamiltonian, 

H = *SJ^(^)\J^) (3.5) 



k 



C " -k I \ C i'-k 



with 



hu> = I M B . n ) , Aa> = 6u>(l - jtfio), B iV = ^ cos k5ii> + \{8 Vji+ i + 5 v ,i-i) (3.6) 



We perform a Bogoliubov transformation by writing, 



k 



b; 
c i,-k 



^(0 +A (O/3|A, fc + r A (O/5 T V fc ) (3-7) 



A>0 



where the /3's obey canonical commutation relations and the two component vectors (f) x (i) = 
(u x (i), v x (i)) are eigenstates of r 3 h, 

r 3 h<j) +x = A0 +A (3.8) 

T 3 h(f)- X = -\(f)- x (3.9) 

Explicitly, <fi~ x = r 1 +A . We normalize the 0's as, 

(0 + V|0 +V > = 5 A)V (3.10) 

Then, up to a constant, 

H = 4 ^EE a ^tVt^ +^ViA, fc ) (3.ii) 

fc A>0 



The solutions to the eigenvalue problem (3.8) with positive eigenvalues can be divided into 



the normalizable and non-normalizable branches. The normalizable branch has dispersion 

A = —7=| sin k y \ (3.12) 
V2 

The continuum branch can be parameterized by momentum < k x < n — k y and has 
dispersion, 

A = \jl — ^(cos k x + cos ky) 2 (3.13) 
We normalize our continuum solutions to, 

(<f>(k x )\T 3 \<f>(k' x )} = (27r)5(k x - k' x ) (3.14) 
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Explicit forms of the eigenstates are given in Appendix [Aj We note that for fixed k y — > 0, 
the energies of both the normalizable state and the continuum threshold tend to ^(A^j, 
with the splitting between these two energies of order k y . This is the reason why the bound 
state does not show up in the effective low energy 0(3) description - it is treated as being 
part of the continuum. 

Now, we can compute the observables. The staggered magnetization is given by, 



TO = S- (c] Cj ) = S- f 72 ^ J2 l^')| 2 ( 3 - 15 ) 



We have evaluated the sum (integral) over the eigenstates numerically - the result is plotted 
in Fig. [3j The staggered moment is depleted near the edge and approaches its bulk value 
monotonically. If we plug 5 = 1/2 into our expansion, the staggered moment at the edge 
is A/edge = 0.217 compared to Nb = 0.303 in the bulk. As already noted, the long distance 
asymptotics of the staggered moment are in good agreement with the predictions of the 0(3) 
continuum theory. 

Similarly, we can compute the bond energies, 

(S 3 S J+X ) = -S 2 + S{(b)b 3 ) + (c] +x c J+x ) + (b 3 c J+x ) + (&t c t +x » 
(S^+y) = -S 2 + S((b]bj) + (c) +y c J+y ) + (b jCj+y ) + (b]c] +y )) 

'^ 2 dk 



S* + S J2mv\j)\ 2 + (u X (j)*v x (j) + v x (j)*u x (j)) cos k y ) 



it/2 K 



A>0 



(3.16) 



The short distance behaviour of the bond energies is shown in Fig. [5j We see that both the 
perpendicular and parallel bonds touching the edge are stronger than in the bulk ((SiSj) is 
more negative), while all the subsequent bonds are weaker than in the bulk. Substituting 
S = 1/2 into our expansion, we find that at the edge {SjSj +x ) = —0.352, (SjSj +y ) = —0.368, 
while in the bulk, (SjSj +fJ ) = —0.329. Thus, comparing to the results of quantum Monte 
Carlo, the 1/S expansion reproduces qualitatively the behaviour of the first two rows of 
bonds away from the edge, but fails to capture the subsequent oscillations in bond strengths 
on short distances. We expect that these oscillations cannot be seen in the perturbative 1/S 
expansion. In the next section, we will argue that the appearance of such oscillations can be 
linked to the existence of a competing valence-bond-solid order parameter. As for the long 
distance asymptotics, we can compare the sum of bond strengths along x and y directions 
to the local energy density computed in the continuum 0(3) model; the two are in good 
agreement (see Fig. Ill) . 
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8<SjS j+fi >/S 



-0.02 



-0.04 



-0.06 



-0.08* 

FIG. 5: Bond strength deviation from bulk value along x (circle) and y (square) directions computed 
in the 1/S expansion. 



Now we turn our attention to the local transverse magnetic susceptibility 

X±Ux) = tJtvt timJ2(S + (jx,qy)S~(j' x ,-q y )) (3.17) 
21 JMy g y ^o 

3x 

where 

S + (j„ q y ) = S 'tl--i»)< ( 3 - 18 ) 

h 

A finite momentum q is needed as a regulator, since we are working in an infinite volume; 
it is convenient to choose q along the y direction. At leading order in the 1/5 expansion, 

X±ti) = ^£<&* + c\_ q ){b], A + c f ,- q )) (3.19) 



= ^ S E E^' 9) + Aj, q))(n X (f, q) + v\f, + 2n(A)) (3.20) 

j> A>0 

where n(X) = (e A//T — 1) _1 is the bose distribution. As expected, for q — > 0, the form-factor 
in (3.20) vanishes upon summing over j', unless A —>■ 0. Thus, we may replace, n(X) — > T/X, 
obtaining, 

x±ti) = s E 9) + vX ^> D)^', q) + « V, ?))* ( 3 - 21 ) 

A short calculation then yields, 



X±U) = ^(1 + (-iy'(V2 + (3.22) 
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This result is saturated by normalizable modes and states at the bottom of the continuum 
band. We see that as j — > oo, the susceptibility approaches its bulk value \x,b — gj- We 
can define the edge susceptibility (per unit edge length) as, 



X_L,edge 



X±,b) 



-3/2 



8J 



(3.23) 



So, at leading order in 1/5* the edge susceptibility is positive, moreover, the approach of x±(j) 
to its bulk value is governed by an oscillating exponential decay. Based on our continuum 
treatment in the previous section, we expect these results to be strongly modified at higher 
orders in 1/S. Indeed, at T = 0, from eq. (2.18) on large distances x±( x ) ~ X±,b falls off as 



1/x. However, the coefficient of this power law is of order 1/5* and, hence, is not captured 



by the leading order result (3.22). When integrated over all space, the large distance power 



law, which is sub leading in the 1/S" expansion, will lead to a logarithmic divergence in the 
size/inverse temperature of the system, which would overpower the leading term in 1/S 



coming from short distances. Thus, the combination of eqs. (2.18), (3.22) naturally explains 



the results of Monte Carlo simulations, which see a positive susceptibility of the "dangling" 
edge spin combined with the negative total edge susceptibility coming from a large distance 
tail in x( x )- 



IV. THE COMB STRUCTURE 

In this section we explain the appearance of the comb structure (Fig. [I]), seen near the 
edge in recent Monte Carlo simulations. In our description, we assume the existence of a 
dynamic valence-bond-solid (VBS) order parameter V(x) with a large correlation length in 
the Neel state. Our treatment becomes exact near a phase transition into a valence-bond- 
solid phase. This phase transition has attracted a lot of attention in the recent years as it lies 
outside the Landau-Ginzburg paradigm.^ It is described by the hedgehog suppressed 0(3) 
cx-model, with the valence-bond- solid order parameter V(x) being the hedgehog insertion 
operator. However, the particular details of the phase transition will not be important for 
our discussion below. 

We begin by defining a microscopic VBS order parameter (which lives on the direct 
lattice), 

V x {i) = (-1)*-+V2 [s{i)S{i + x)- S{i)S{i - x)) (4.1) 
V y {i) = (-l)^ 2 (s(i)S(i + y)- S(^S(z - y)) (4.2) 

In this section, we take the origin to lie on the dual lattice. It is customary to group V x , 
V y into a complex order parameter V — V x + iV y which has the following transformation 
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a) b) 

FIG. 6: a) Lattice order with (V x ) / 0. b) Bond order with (O x ) / 0. 



properties under elements of the square lattice space group: 



y;.n /,../,)•/:, = -v\i x -i } i y ) (4.3) 

TlV{l X ,ly)Ty = V\l X ,ly-l) (4.4) 

I^V{i x ,i y )I^ = V{ _ lx ^ ly) (4.5) 

ll duaX V{i x ,i y )I^ = V(i x ,-i y ) (4.6) 

Rlt Xv ^h)K72 = iV\i y -i x ). (4.7) 



Here T XtV are translations by one lattice spacing in the x, y directions, I xy are x, y, reflections 
about a dual lattice point, and R^/f is a 90° rotation about about a dual lattice point. For 
completeness we also list the transformation property of V under rotations about direct 
lattice point (-1/2,-1/2), 

R[f 2 V{i x , i y )R% = iV(i y , -1 - i x ) (4.8) 

A non-zero expectation value of the VBS order parameter V would lead to a bond pattern 
shown in Fig. [6] a). As already noted, the operator V(x) is represented by the hedgehog 
insertion operator in the continuum description of the antiferromagnet - valence bond solid 
transition. 

Clearly, the order parameter V is adequate for describing the oscillations of horizontal 
bonds in the comb structure (Fig. [T|. However, the oscillations of the vertical lines in the 
comb structure (Fig. [I]), shown separately in Fig. [6] b) are not of the "dimer form." To 
describe them, we introduce a new order parameter, 
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O x {i) = \S{i +\x+ l -y)S(i + ^x- l -y) - S{i - l -x + l -y)S{i - X -x - l -y) 

O y {i) = (-l) iy (s(i + ^y+ l -x)S{i + l -y- ^x) - S(i - l -y+ l -x)S{i - X -y - l -x] 

O x describes vertical bond lines which are oscillating in strength along the x direction (see 
Fig. [6]b)). Similarly, O y describes horizontal bond lines, which are oscillating in strength 
along the y direction. 

We can group O x and O y into a single complex order parameter O = O x + iO y . The 
transformation properties of O are, 

r jrn.,.. i.vr,. = -o\i x -i,i y ) (4.9) 

TlO{l X ,ly)Ty = 0\l X ,ly-l) (4.10) 
lt^O{i x ,iy)I^ = -0\-i x ,iy) (4.11) 
4 dUaX 0(i x ,iy)I^ = 0\i X ,-iy) (4.12) 

Rl% al O(t x ,i y )R^f = iO(i y ,-Q (4.13) 

and for rotations about direct lattice point (—1/2, —1/2): 

Rlf 2 0(z x , i y )R% = iO\i v , -1-4) (4.14) 

Now we may ask whether it is possible in the continuum to construct an operator with the 
transformation properties of 0(x) out of V(x). Clearly, any function of V with no derivatives 
cannot do the job, since under dual lattice reflections I x u y al , O transforms non-trivially, while 

V transforms trivially. Thus, a static uniform condensate of V (not surprisingly) cannot give 
rise to the order in Fig. [6]b). However, we can obtain an expression with the transformation 
properties of O if we allow for derivatives of V. Considering expressions with one power of 

V and one derivative, we obtain, 

X ~ 8 X V X , Oy ~ dyVy (4.15) 

(with the same proportionality constant). 

Thus, if dimerization of horizontal bonds is present and is inhomogeneous along the x 
direction then we automatically obtain the "secondary" order in Fig. |6]b). 

Now, we may ask, how a non-zero expectation value of the VBS order is generated? 
Indeed, in the Neel phase, in the bulk, the Z4 lattice rotation symmetry is unbroken and 
(V) = 0. However, the edge possesses a smaller lattice symmetry group than the bulk - 
in particular, the lattice rotation symmetry is explicitly broken. This is manifested in the 
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continuum formulation by the appearance of an edge perturbation, 



5S = hi J drdy (V + V f ) = h J drdy V x (4.16) 

In the phase where V is gapped, we expect such a coupling will lead to an appearance of 
(V x (x, y)) decaying away from the edge. Hence, we will also have (O x (x, y)) ^ 0, which close 



to the critical point can just be obtained from (4.15). Thus, the appearance of the comb 
structure is very natural. 

Based on the known results on boundary critical behaviour, 1 ^ we may write down the 
scaling forms for (V(x)), (0(x)) in the critical region. The edge perturbation 5S is relevant 
at the critical point provided that A v < 2, where A v is the scaling dimension of operator 
V(x). Then the scaling forms become universal (up to overall multiplicative factors), 

<Y*{x)) ~ ^U(*/0 (4-17) 
(O x (x)) ~ ^g'ix/O (4.18) 

Here £ is the correlation length of the VBS order parameter in the Neel phase (which 
is proportional to the inverse spin stiffness c/p s with some universal amplitude). In the 
deconfined criticality scenario, £ will be given by the inverse skyrmion mass. Note that 
due to the extra derivative in O compared to V, the modulations of lines parallel to the 
edge become parametrically weaker than those of dimers perpendicular to the edge as we 
approach the phase transition. We may also write down short and long distance asymptotics 
of g(u), 

g{u) ~ ^v, u^O (4.19) 



u 



g{u) ~ e~ u , u — > oo (4.20) 

where we have not specified the likely power-law prefactor for the long distance asymptotic 
P~20l). 



V. A NON-MAGNETIC IMPURITY 

This section will briefly discuss the case of a different defect in a perfect square lattice 
antiferromagnet: a single site with a missing spin. This is often experimentally realized in 
Cu antiferromagnets by replacing Cu with Zn. We are interested in the configuration of VBS 
order around this impurity — this was addressed recently in Ref. |9] using methods similar to 



those used in Section IV Our purpose is to connect these phenomenological approaches to 



the field-theoretic treatment near the deconfined critical point presented in Ref. [TUJ. 



As in Section IV, we begin by describing the influence of the impurity by writing down the 
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action for V consistent with the symmetries of the impurity Hamiltonian; here the action 
has to be invariant under i?^ *? -^ irec \ and time-reversal. Then the analog of the edge 



perturbation in Eq. (4.16) for an impurity at Xi mp 1 



imp, V 



-Ai / dr 



dV dV^ M _ 

dx dx dy dy 



(5.1) 



As shown in Ref . EJ this perturbation induces 'vortices' in the VBS order around the impurity. 
We now discuss the origins of the term <Si mPi y in the critical theory of the Neel-VBS transition 
in the insulator. This will determine the behavior of the coupling Ai near this transition. 

The behavior of a non-magnetic impurity near this transition has been described in 
Refs. [T0)13fl4l For the bulk model without an impurity, a field theoretic description of the 



vicinity of the quantum critical pointP^^^ is provided by the 



N-l 



theory at N = 2: 



S 



d rdr 



- iA^)z a \ 2 + s\z a \ 2 + | (\z a \ 2 ) 2 + 



2e 2 



d u A x f 



(5.2) 



Here /x, v, A are spacetime indices, z a , a = 1 . . . N = 2 is a complex scalar which is a 
SU(iV) fundamental, and is a non-compact U(l) gauge field. As discussed in Ref. 



the most important perturbation to Eq. (5.2) induced by the non-magnetic impurity near 



the deconfined critical point is the impurity Berry phase: 



5 



imp 



iQ / drA T {x = 0,r) 



(5.3) 



where Q is a 'charge' characterizing the impurity. The value of Q does not flow under the 
RG, and so Q is a pure number which controls all universal characteristics of the impurity 
response. 

Let us now discuss the symmetries of iS+iSi mp . In addition to the global SU(iV) symmetry, 
this model has a global U(l) v symmetry which is the dual of the U(l) gauge invariance. The 
primary action of this symmetry is on the monopole operator, V(x,r), which transforms as 



U(l)„ : V 
Ro : V 



Ve lLp 
V 



(5.4) 



At the moment, this U(l)^ 'flux' symmetry is independent of spatial rotations Rg, and this 
has been indicated above for completeness. The physical Z 4 lattice rotation symmetry is the 
combination of tt/2 rotations in U(l) v and Rg - thus, the monopole operator V is identified 
with the VBS order parameter. 

A key property^ of the theory <S+iSi mp is the operator product expansion for the monopole 
operator V in the vicinity of the impurity 



lim V(x, r 
l*l—o 



~ \x\ im p e 



^imp(^~) 



(5.5) 
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where 9 is the azimuthal angle of x, and AjjL, is the impurity correction to the scaling 



imp 



dimension of V (A y ) as defined in Ref. [TUl Here VJ mp is a fluctuating impurity degree of 
freedom with a non-trivial scaling dimension. The presence of the e~ l ® d factor indicates a 
Q-fold winding in the phase of the VBS order parameter around the impurity. Thus, the 
effect of the impurity Berry phase term is to induce vortex-like correlations in bond order 
near the impurity. 

However, the way this vortex is pinned to the lattice is determined by additional impurity 



perturbations, the most relevant of which is given by Eq. (5.1). We can understand this 
by continuing our symmetry analysis. The combination of Eqs. (5.4) and (5.5) implies the 
following transformations of Vj mp under the flux symmetry and spatial rotations: 

Tjm ■ y. — > V- e ilfi 

yJ \- L Jtp • 'imp 7 'imp 

Re : V imp -> V imp e- iB (5.6) 

Here, and henceforth, we specialize to the case Q — 1, although the generalization to other 
Q is not difficult. We note that the quantum numbers of Vj mp are the same as those of the 



perturbation (5.1). Hence, the two will mix and we may replace (5.1) by, 



SL P ,v = -Kj drV imp (r) (5.7) 



Now, there are two possibilities. If the perturbation (5.7) is relevant at the critical point, 
which occurs for 

dim[M mp ] = A V + A£ p < 1 (5.8) 

the coupling X[ will flow to infinity. In this case, at criticality, the VBS order parameter will 
be given by, 



(y&T))~J^y (5.9) 

\x A 



Alternatively, if the coupling X[ is irrelevant, we can treat it in perturbation theory and 
obtain, 

( V &r))~ j" A y ■ (5.10) 

\0C\ imp— 1 

Now let us move into the Coulomb phase of S, where there is a mass gap, m, for the z a 
spinons. We are interested in the effective theory for V(x) at energy scales smaller than this 
mass gap. The only low energy degree of freedom is the (pseudo)-Goldstone (p associated 
with spontaneous breaking of the 11(1)^ symmetry. We identify, V ~ m AV e lv . The effective 
action for the tp field in the absence of impurity takes the form, 

S = J d 2 xdr (^-^(d^f - A 4 cos(4<^ (5.11) 
Here e ~ m is the effective electric charge in the Coulomb phase and A4 ~ m dim[y4] is the 
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dangerously irrelevant perturbation that breaks U(l) v symmetry to the physical Z4. Now, 
let's discuss the impurity perturbations in the effective theory. One such perturbation can 
be simply obtained from (5.1) by replacing V — > e % 



Ss 



imp,efT 



— Ai iG ff J dr (i(d x ip + id y (p)e %Lp + h.c.) 



(5.12) 



We are interested in how the coefficient of this term Ai i6 g is renormalized. If the perturbation 



(|5.ip is relevant (see eq. ([578])) , the impurity response will be universal and Xi >e s will be a 

1 



constant in m by dimensional analysis. Otherwise, Ai )e fj 



A V +AY - 
~ 171 lm P 



VI. CONCLUSION 

In this paper we have addressed two puzzles raised by recent Monte Carlo simulations of 
edge response in square lattice quantum antiferromagnets. The first puzzle is the appearance 
of negative edge susceptibility - we have shown that this effect is due to low energy spin- 
waves. We predicted that the total edge susceptibility diverges logarithmically as inverse 
temperature/system size goes to infinity, and found this to be in good agreement with the 
Monte Carlo simulations of Ref. El We would like to note here that our results on the low 
temperature behaviour of susceptibility apply equally well to a clean and rough edge, as 
our continuum 0(3) cr-model description does not assume translational invariance along the 
edge. (However, for the rough edge, there may be additional important contributions to the 
susceptibility coming from Berry phase effects, not present in the 0(3) a model.) The second 
puzzle is the observation of a comb structure in the bond response near the edge. We have 
argued that this is likely a purely quantum mechanical effect, which cannot be captured 
by the naive 1/S expansion. We have shown that the appearance of the comb structure 
can be understood in the framework of a continuum theory involving a dynamical valence- 
bond-solid order parameter. Such a description becomes exact in the neighbourhood of a 
quantum phase transition to a valence-bond-solid phase. We hope that the simulations of 
edge response in Heisenberg modeP will be extended to the so-called JQ model where such 
a phase transition is observed.^ We have made a few predictions regarding the behaviour of 
the comb structure near criticality, e.g. the relation between the behaviour of bonds parallel 
and perpendicular to the edge in the comb. Edge response near the quantum critical point 
might also be a viable way to extract the scaling dimension of the valence-bond-solid order 



parameter, see eqs. (4. 17), (4. 19) 



Finally, in Section |Vj we briefly discussed some related issues on the problem of a single 
non-magnetic impurity, complementary to the more detailed discussion of this case in Ref. El 
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APPENDIX A: EIGENFUNCTIONS OF BOGOLIUBOV QUASIPARTICLES 

First, we define for fixed energy A, 

1 / -v / TTA 



u 

V 



(Al) 



" V2X \ VI- A 

Now, the eigenstates can be expressed as, 

Normalizable solution: 

A = -^=\smk y \, -7r/2 < k y < rr/2 (A2) 
V2 

<f>(j) = Cl j + c 2 ^ ~ M j (-lYe-™ (A3) 
e s i = (y/2 + l)( v / l + cos 2 A; y - cos k y ), e S2 = (V2 + l)(^/l + cos 2 fcj, + cos fc„) 



Cl 



2-t(v^-l)|sinA; w | / e S2 V / T^A - V / TTA 



c 2 i - | sin ky\y/l + cos 2 fcj, I e 8l y/l -X-y/l + X 



(A4) 



Continuum solutions: 



7 = ^(cos(A; x ) + cos(fcj,)), A = a/1 — 7 2 , < < 7r — \k y \, —ir/2 < k y < n/2 

(A5) 

Branch 1: < k x < cos _1 (l — 2cos/c J/ ) 

1 ( (ae ikxj + a*e~ ikxj — (—l) j e~ js )u \ ^ 
<t>{j) = o « , * -« 7 , 1 7 -7, > s = cosh- 1 (cosA; a; + 2cosA; ?/ ) (A6) 
\a\ \ (ae %kxl + a*e lkx3 + (-lye JS )v I 



a = I 7e s - 1 - -— -((jcosk x - l)e s + 7 - cosk x ) 

2A V sin/^ 



(A7) 
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Branch 2: cos x (l — 2cosk y ) < k x < 7T — \k y \ 



k x = 7r — cos (cos(k x ) + 2 cos(ky)) , ir — \k y \ < k x < n (A£ 
o ..{./) - -4 | r- LJ (■(»«( /,',(./- 1/21) ( " ) H-CiaCOS^O' + 1/2)) ( ~ U I ) (A9) 




2 (j) = A ( c 21 sin(Ai(j + I/ 2 )) ( " J + c 22 sin(fc,(j + 1/2)) ( " ) ) (A10) 

A= (sinfc ;E )t(sin((fc :c + fc :c )/2)+7sin((A; ;c -A; ;c )/2))-5 (All) 

,i f2coskJ2\* . .i f2cosk x /2\^ /A , 

cii = (1 + 7 M t4H , cu = (l-7 M A12 

y cosfc x /2 y V cos/^/2 / 

.'2sin^/2\ 3 . ,i /2sinA; ;c /2\ ^ 

cai = 1-7 ' ■ ,1 , Qa = (l + 7 M . r a ; o A13 
sin^/2 y V sinfc^/2 / 

The division of the continuum spectrum into two branches is clear when we look at a plot 



of X(k x ) (A5) for k y fixed. For Wl — cos 4 (k y /2) < X < y 1 — sin (k y /2) there is only one 



corresponding value of k x in the range < k x < ir (there is always a solution with opposite 



k x , as well). This is our branch 1. On the other hand, for yl — sin 4 (fc y /2) < A < 1 there 
are two solutions with < k x < ir, which we label by k x and k x . These two solutions are 
mixed by the edge and form the two linearly independent eigenstates <fii, <p2 in branch 2. 
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